Method for enhancing blood vessels in angiography images

ABSTRACT

Disclosed is a method for enhancing blood vessels in angiography images. The method incorporates the use of linear directional features present in an image, extracted by a Directional Filter Bank, to obtain more precise Hessian analysis in noisy environment and thus can correctly reveal small and thin vessels. Also, the directional image decomposition helps to avoid junction suppression, which in turn, yields continuous vessel tree.

BACKGROUND

1. Technical Field

The invention generally concerns enhancement filtering to improvevisibility of blood vessels and more practically to a framework forvessel enhancement filtering in angiography images.

2. Description of the Related Art

The common way to interpret vasculature images, e.g. the MagneticResonance Angiography (MRA) images, is to display them in MaximumIntensity Projection (MIP) in which the stack of slices is collapsedinto a single image for viewing. MIP is performed by assigning to eachpixel in the projection the brightest pixel over all slices in thestack. With this type of display, small vessels with low contrast arehardly visible and other organs may be projected over the arteries. FIG.1 may demonstrate that small vessels tend to resemble background. Avessel enhancement procedure as a pre-processing step for maximumintensity projection display will help to diminish these twolimitations.

There are a variety of vessel enhancement methods in literature. Thesimplest one is to threshold the raw data but this makes thesegmentation process incorrectly label bright noise regions as vesselsand cannot recover small vessels which may not appear connected in theimage. Recently, Hessian-based approaches have been utilized in numerousvessel enhancement filters. These filters are based on the principalcurvatures, which are determined by the Hessian eigenvalues, todifferentiate the line-like (vessel) from the blob-like (background)structures. However, their disadvantage is that they are highlysensitive to noise due to second-order derivatives. Moreover, they tendto suppress junctions which are characterized same as the blob-likestructures using the principal curvature analysis. Junction suppressionin turn leads to discontinuity of the vessel network.

SUMMARY

The present invention has been made to solve the above problemsoccurring in the prior art. There is provided a method for the vesselenhancement filter utilizing the linear directional information presentin an image. The method comprises decomposing the input angiographyimage into directional images T_(i) using Decimation-free DirectionalFilter Bank (DDFB), removing non-uniform illumination by employing ndistinct homomorphic filters matched with its corresponding directionalimage, enhancing vessels in every directional image, and re-combiningall enhanced directional images. Further consistent with the presentinvention, wherein said DDFB comprises filtering the input angiographyimage with H₀₀(ω₁, ω₂) and H₁₁(ω₁, ω₂) hourglass-shaped like passbands,filtering with H₀₀(Q^(T)(ω₁, ω₂)) and H₁₁(Q^(T)(ω₁, ω₂)), where Trepresents transpose and Q is Quincunx downsampling matrix, andfiltering with H₀₀(R_(i)Q^(T)Q^(T)(ω₁, ω₂)) and H₁₁(R_(i)Q^(T)Q^(T)(ω₁,ω₂)) where R_(i) (i=1, 2, 3, and 4) are resampling matrices.

$Q = \begin{pmatrix}1 & 1 \\{- 1} & 1\end{pmatrix}$ $R_{1} = \begin{pmatrix}1 & 1 \\0 & 1\end{pmatrix}$ $R_{2} = \begin{pmatrix}1 & {- 1} \\0 & 1\end{pmatrix}$ $R_{3} = \begin{pmatrix}1 & 0 \\1 & 1\end{pmatrix}$ $R_{4} = \begin{pmatrix}1 & 0 \\{- 1} & 1\end{pmatrix}$

Output of the vessel enhancement filter for one directional image is

${{\Phi (p)} = {\max\limits_{\sigma \in S}{\varphi_{\sigma}(p)}}},$

where p is coordinate (x′,y′), S is a range, and σ is a various scale.The coordinates Ox′y′ is obtained by rotating Oxy by the angleassociated with that directional image. φ_(σ)(p) is based on thediagonal values of the Hessian matrix H′ in the coordinates Ox′y′.

$H^{\prime} = \begin{bmatrix}\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} & \frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} \\\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} & \frac{\partial^{2}I_{i}}{\partial y^{\prime 2}}\end{bmatrix}$ where${\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\cos^{2}\theta_{i}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\sin^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{\partial y^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\cos^{2}\theta_{i}} - {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\cos^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} = {{{- \frac{1}{2}}\frac{\partial^{2}I_{i}}{\partial x^{2}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\cos \left( {2\theta_{i}} \right)}} + {\frac{1}{2}\frac{\partial^{2}I_{i}}{\partial y^{2}}{\sin \left( {2\theta_{i}} \right)}}}}$

Specifically, the input image is first decomposed by DDFB into a set ofdirectional images, each of which contains linear features in a narrowdirectional range. The directional decomposition has two advantages. Oneis, noise in each directional image will be significantly reducedcompared to that in the original one due to its omni-directional nature.The other is, because one directional image contains only vessels withsimilar directions, the principal curvature calculation in it isfacilitated. Then, distinct appropriate enhancement filters are appliedto enhance vessels in the respective directional images. Finally, theenhanced directional images are re-combined to generate the output imagewith enhanced vessels and suppressed noise. Thisdecomposition-filtering-recombination scheme also helps to preservejunctions.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is showing an angiography image with small vessels.

FIG. 2 is flowchart showing a method of enhancing blood vesselsconsistent with the present invention.

FIGS. 3A and 3B are showing the frequency responses of hourglass-shapedlike filters.

FIG. 4A is showing the First stage of DDFB structure.

FIG. 4B is showing the Second stage of DDFB structure.

FIG. 4C is showing the Third stage of DDFB structure.

FIG. 5 is showing the block diagram of the present invention.

FIGS. 6A, 6B, 6C, 6D, 6E, 6F, 6G, and 6H are showing eight demonstrativedirectional images.

FIG. 7A is showing a synthetic image used to evaluating the performanceof the present invention.

FIGS. 7B, 7C, and 7D are showing enhancement results of the Frangifilter, the Shikata filter and the present invention for the input imageshown in FIG. 7A.

FIGS. 8A, 8B, and 8C are showing enhancement results of the Frangifilter, the Shikata filter and the present invention for the input imageshown in FIG. 1.

DETAILED DESCRIPTION Embodiment 1

Whereinafter, a embodiment consistent with the present invention will bedescribed with reference to the drawing.

The proposed method consists of three steps, as shown in FIG. 2: Firststep (step 21) is construction of directional images by decomposinginput image, second step (step 22) is vessel enhancement, and third step(step 23) is recombination of enhanced directional images.

As shown in FIG. 2, the decomposing with DDFB (step 21), which is thefirst step of decomposing input angiography image into directionalimages. Next, removing non-uniform illumination by homomorphic filterand enhancing directional images by appropriate enhancement filters(step 22), which is the second step of enhancement filtering to improvevisibility of blood vessels. Thereafter, re-combining directional images(step 23), which is the third step of re-combining all enhanceddirectional images.

So, the invention is characterized as follows.

To enhance vessels in angiography images, input vessel image isdecomposed to a set of directional images using DDFB. The non-uniformillumination is removed by employing a homomorphic filter matched withits corresponding directional image. The filtering process is based onthe Hessian eigenvalues and filtering process is applied on the set ofdirectional images.

Directional Filter Bank (DFB) can decompose the spectral region of aninput image into n=2^(k) (k=1, 2, . . . ) wedge-shaped like subbandswhich correspond to linear features in a specific direction in spatialdomain.

One disadvantage of DFB is that the subbands are smaller in size ascompare to the size of input image. The reduction in size is due to thepresence of decimators. As far as image compression is concerned,decimation is a must condition. However, when DFB is employed for imageanalysis purposes, decimation causes two problems. One is, as weincrease the directional resolution, spatial resolution starts todecrease, due to which we loose the correspondence among the pixels ofDFB outputs. The other is, an extra process of interpolation is involvedprior to enhancement implementation. This extra interpolation processnot only affects the efficiency of whole system but also produces falseartifacts which can be harmful especially in case of medical imagery.Some vessels may be broken and some can be falsely connected to someother vessels due to the artifacts produced by interpolation. So a needarises to modify directional filter bank structure in a sense that nodecimation is required during analysis section. We suggest to shift thedecimators and resamplers to the right of the filters to create theDDFB, which yields directional images rather than directional subbands.This consequently results in elimination of interpolation and naturallyfits the purposes of feature analysis.

The decomposing step (step 21) of applying DDFB comprise as followingstages. First stage of filtering the input angiography image withH₀₀(ω₁, ω₂) and H₁₁(ω₁, ω₂) hourglass-shaped like passbands, Secondstage of filtering with H₀₀(Q^(T)(ω₁, ω₂)) and H₁₁(Q^(T)(ω₁, ω₂)), whereT represents transpose and Q is Quincunx downsampling matrix, and Thirdstage of filtering with H₀₀(R_(i)Q^(T)Q^(T)(ω₁, ω₂)) andH₁₁(R_(i)Q^(T)Q^(T)(ω₁, ω₂)).

At first the stage of applying DDFB, construction of first stage of DDFBonly requires two filters. Filters at first stage of DDFB are H₀₀(ω₁,ω₂) and H₁₁(ω₁, ω₂). They have hourglass-shaped like passbands as shownin FIGS. 3A and 3B. FIG. 4A shows the block diagram of the first stageof DDFB.

At second the stage of applying DDFB, the filters required forconstruction of second stage are H₀₀(Q^(T)(ω₁, ω₂)) and H₁₁(Q^(T)(ω₁,ω₂)), where T represents transpose and Q is the Quincunx downsamplingmatrix.

$\begin{matrix}{Q = \begin{pmatrix}1 & 1 \\{- 1} & 1\end{pmatrix}} & {{EQUATION}\mspace{20mu} 1}\end{matrix}$

Spectral regions of directional images obtained after filtering throughsecond stage filter are shown in FIG. 4B.

At third the stage of applying DDFB, filters used during the third stageof DDFB are H₀₀(R_(i)Q^(T)Q^(T)(ω₁, ω₂)) and H₁₁(R_(i)Q^(T)Q^(T)(ω₁,ω₂)), as shown in FIG. 4C where R_(i) (i=1, 2, 3, and 4) are resamplingmatrices.

$\begin{matrix}{{R_{1} = \begin{pmatrix}1 & 1 \\0 & 1\end{pmatrix}}{R_{2} = \begin{pmatrix}1 & {- 1} \\0 & 1\end{pmatrix}}{R_{3} = \begin{pmatrix}1 & 0 \\1 & 1\end{pmatrix}}{R_{4} = \begin{pmatrix}1 & 0 \\{- 1} & 1\end{pmatrix}}} & {\text{-}{EQUATION}\mspace{20mu} 2\text{-}}\end{matrix}$

Overall eight different filters are created to be used during the thirdstage.

By using the DDFB, the input image is decomposed to n=2^(k) (k=1, 2, . .. ) directional images T_(i). The motivation here is to detect thin andlow-contrast vessels (which are largely directional in nature) whileavoiding false detection of non-vascular structures. Directionalsegregation property of DDFB is helpful in eliminating randomly orientednoise patterns and non-vascular structures which normally yield similaramplitudes in all directional images (see FIGS. 6A to 6H).

Before these directional images are enhanced in the next step, they areutilized to efficiently remove non-uniform illumination (NUI), whichlimits the correct vessel enhancement. One conventional approach tocorrect NUI is to directly apply homomorphic filtering on the originalimage. A general image can be characterized by two components: (1) theillumination component, which changes slowly in a neighborhood due tolight source characteristics and thus is low-frequency, and (2) thereflectance component, which is determined by the amount of lightreflected by the objects and therefore is high-frequency. Thehomomorphic filter is to suppress the low-frequency component whileenhance the high-frequency one. However, the direct application ofhomomorphic filtering is sometimes unsatisfactory because it may enhancebackground noise which is normally high-frequency. Tuning the filterparameters in this case suffers from a compromise. The more NUI isremoved, the more background noise is enhanced. Differently, we proposeemploying a homomorphic filter matched with its correspondingdirectional image as shown in the dash-boundary box in FIG. 5. This newarrangement provides us a better control on the parameters of individualhomomorphic filter.

Explaining the second step (step 22) of vessel enhancement, we proposeremoving non-uniform illumination by homomorphic filter.

In order to compute the principal curvatures with less noisesensitiveness, it is necessary to align the vessel direction with thex-axis. One possible way is to rotate the directional images.Nevertheless, image rotation requires interpolation which is likely tocreate artifacts and thus is harmful especially in case of medicalimagery. We therefore rotate the coordinates rather than the directionalimages.

Suppose the directional image I_(i) (i=1 . . . n) corresponds to theorientations ranging from θ_(i,min) to θ_(i,max) (counterclockwiseangle). Its associated coordinates Oxy will be rotated to Ox′y′ by anamount as large as the mean value θ_(i).

$\begin{matrix}{\theta_{i} = \frac{\theta_{i,\min} + \theta_{i,\max}}{2}} & {\text{-}{EQUATION}\mspace{20mu} 3\text{-}}\end{matrix}$

Hessian matrix of the directional image I_(i) in the new coordinatesOx′y′ is determined as followed EQUATION 4.

$\begin{matrix}{{H^{\prime} = \begin{bmatrix}\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} & \frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} \\\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} & \frac{\partial^{2}I_{i}}{\partial y^{\prime 2}}\end{bmatrix}}{where}{\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\cos^{2}\theta_{i}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\sin^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{\partial y^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\sin^{2}\theta_{i}} - {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\cos^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} = {{{- \frac{1}{2}}\frac{\partial^{2}I_{i}}{\partial x^{2}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\cos \left( {2\theta_{i}} \right)}} + {\frac{1}{2}\frac{\partial^{2}I_{i}}{\partial y^{2}}{\sin \left( {2\theta_{i}} \right)}}}}} & {\text{-}{EQUATION}\mspace{20mu} 4\text{-}}\end{matrix}$

The principal curvatures are then defined by the diagonal values of H′.These values are EQUATION 5.

$\begin{matrix}{{{{PC}_{1} = 0};}{{PC}_{2} = {\frac{y^{\prime 2} - \left( {\sigma_{0}^{2} + \sigma^{2}} \right)}{\left( {\sigma_{0}^{2} + \sigma^{2}} \right)^{2}}{I_{i}\left( {x^{\prime},y^{\prime}} \right)}}}} & {\text{-}{EQUATION}\mspace{20mu} 5\text{-}}\end{matrix}$

where σ selected in a range S is the standard deviation of the Gaussiankernel used in the multiscale analysis.

Practically, the vessel axis is not, in general, identical to thex′-axis and so PC₁≈0.

Inside the vessel, |y′|<√{square root over (σ₀ ²+σ²)} and thus PC₂ isnegative. Therefore, vessel pixels are declared when PC₂<0 and

${\frac{{PC}_{1}}{{PC}_{2}}}\mspace{11mu} \text{<<}\mspace{11mu} 1.$

To distinguish background pixels from others, we define a structurenessmeasurement as EQUATION 6.

C=√{square root over (PC ₁ ² +PC ₂ ²)}  EQUATION 6

This structureness C should be low for background which has no structureand small derivative magnitude.

Based on the above observations, the vessel filter output can be definedas EQUATION 7.

$\begin{matrix}{{{\varphi_{\sigma}(p)} = {{\eta \left( {PC}_{2} \right)}{{\exp\left( {- \frac{R^{2}}{2\beta^{2}}} \right)}\left\lbrack {1 - {\exp\left( {- \frac{C^{2}}{2\gamma^{2}}} \right)}} \right\rbrack}}},} & {\text{-}{EQUATION}\mspace{20mu} 7\text{-}}\end{matrix}$

where p=(x′,y′), R=PC₁/PC₂, β and γ are adjusting constants, and

${\eta (z)} = \left\{ \begin{matrix}0 & {{{{if}\mspace{14mu} z} \geq 0};} \\1 & {{{if}\mspace{14mu} z} < 0.}\end{matrix} \right.$

The filter is analyzed at different scales σ in a range S. When thescale matches the size of the vessel, the filter response will bemaximum. Therefore, the final vessel filter response is EQUATION 8.

$\begin{matrix}{{\Phi (p)} = {\max\limits_{\sigma \in S}{\varphi_{\sigma}(p)}}} & {\text{-}{EQUATION}\mspace{20mu} 8\text{-}}\end{matrix}$

One filter (EQUATION 8) is applied to one directional image to enhancevessel structures in it.

Explaining the third step (step 23) of re-combining directional images,each directional image now contains enhanced vessels in its directionalrange and is called the enhanced directional image.

Denote Φ_(i)(p), i=1 . . . n, as the enhanced directional images.Another advantage of DDFB is that its synthesis is achieved by simplysumming all directional images. Thus, the output enhanced image F(p) canbe obtained by EQUATION 9.

$\begin{matrix}{{F(p)} = {\frac{1}{n}{\sum\limits_{i = 1}^{n}{\Phi_{i}(p)}}}} & {\text{-}{EQUATION}\mspace{20mu} 9\text{-}}\end{matrix}$

The whole filtering procedures can be summarized as follows. First, theinput angiography image is decomposed into n=2^(k) (k=1, 2, . . . )directional images T_(i) using DDFB. Then, n distinct homomorphicfilters are employed to n respective directional images to removenon-uniform illumination. The output uniformly illuminated directionalimages I_(i) are enhanced according to EQUATION 7 and EQUATION 8.Finally, all enhanced directional images are re-combined to yield thefinal filtered image F as in EQUATION 9.

FIG. 7 shows the results of an synthetic image which was processed bythe three filter models. The synthetic image in FIG. 7A is designed tocontain vessels of different sizes and junctions of different types. Itis possible to see that the Frangi (FIG. 7B) and Shikata (FIG. 7C)filters unexpectedly suppress junctions while our proposed approach(FIG. 7D) does not. The suppressed junctions make vessels discontinuous.

It is the use of directional image decomposition that makes the proposedmodel work. Normally, a vessel has one principal direction, which ismathematically indicated by a small ratio between the smaller and largerHessian eigenvalue. Meanwhile, at a junction, where a vessel branchesoff, there are more than two principal directions, and thus the ratio oftwo eigenvalues is no longer small. As a result, the conventionalenhancement filters (e.g., the Frangi and Shikata filters) considerthose points as noise and then suppress them. Our proposed approach, onthe other hand, decomposes the input image to various directionalimages, each of which contains vessels with similar orientations.Consequently, junctions do not exist in directional images and so arenot suppressed during the filtering process. After that, due to there-combination of enhanced directional images, junctions arere-constructed at those points which have vessel values in more than twodirectional images.

FIGS. 8A, 8B, and 8C respectively show enhancement results of Frangifilter, Shikata filter and our present invention for the input imagesshown in FIG. 1. As can be observed, Frangi filter gives good resultswith large vessels but fails to detect small ones while Shikata model isable to enhance small vessels but unfortunately enhances backgroundnoise also. Conversely, our proposed filter can enhance small vesselswith more continuous appearances.

1. A method for enhancing blood vessels in angiography image, the methodcomprising the steps of: decomposing the angiography image into ndirectional images T_(i) using DDFB; removing non-uniform illuminationby employing n distinct homomorphic filters matched with itscorresponding directional image; enhancing n directional images using nvessel enhancement filters; and re-combining all enhanced directionalimages.
 2. The method as claimed in claim 1, wherein the number n ofdirectional images T is n=2^(k) (k=1, 2, . . . ).
 3. The method asclaimed in claim 1, wherein said decomposing the angiography imagecomprises steps of; filtering the input angiography image with H₀₀(ω₁,ω₂) and H₁₁(ω₁, ω₂), each having hourglass-shaped like passbands;filtering with H₀₀(Q^(T)(ω₁, ω₂)) and H₁₁(Q^(T)(ω₁, ω₂)), where Trepresents transpose and Q is Quincunx downsampling matrix; andfiltering with H₀₀(R_(i)Q^(T)Q^(T)(ω₁, ω₂)) and H₁₁(R_(i)Q^(T)Q^(T)(ω₁,ω₂)), where i=1, 2, 3, 4, and R_(i) are $R_{1} = \begin{pmatrix}1 & 1 \\0 & 1\end{pmatrix}$ $R_{2} = \begin{pmatrix}1 & {- 1} \\0 & 1\end{pmatrix}$ $R_{3} = \begin{pmatrix}1 & 0 \\1 & 1\end{pmatrix}$ $R_{4} = {\begin{pmatrix}1 & 0 \\{- 1} & 1\end{pmatrix}.}$
 4. The method as claimed in claim 3, wherein R_(i) andQ are respectively resampling and downsampling matrices.
 5. The methodas claimed in claim 1, wherein output of the vessel enhancement filteris${{\Phi (p)} = {\max\limits_{\sigma \in S}{\varphi_{\sigma}(p)}}},$where p is coordinate (x′,y′), S is a range, and σ is various scale. 6.The method as claimed in claim 5, wherein the coordinates Ox′y′ isobtained by rotating Oxy by an angle defined as mean value θ_(i), where$\theta_{i} = \frac{\theta_{i,\min} + \theta_{i,\max}}{2}$
 7. The methodas claimed in claim 5, wherein the vessel enhancement filer at a certainscale σ is${{\varphi_{\sigma}(p)} = {{\eta \left( {PC}_{2} \right)}{{\exp\left( {- \frac{R^{2}}{2\beta^{2}}} \right)}\left\lbrack {1 - {\exp\left( {- \frac{C^{2}}{2\gamma^{2}}} \right)}} \right\rbrack}}},$where PC₁ and PC₂ are the diagonal values of the Hessian matrix H′ inthe coordinates Ox′y′, R=PC₁/PC₂, C=√{square root over (PC₁ ²+PC₂ ²)}, βand γ are adjusting constants, and ${\eta (z)} = \left\{ \begin{matrix}0 & {{{{if}\mspace{14mu} z} \geq 0};} \\1 & {{{if}\mspace{14mu} z} < 0.}\end{matrix} \right.$
 8. The method as claimed in claim 7, wherein theHessian matrix is ${H^{\prime} = \begin{bmatrix}\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} & \frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} \\\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} & \frac{\partial^{2}I_{i}}{\partial y^{\prime 2}}\end{bmatrix}},{where}$${\frac{\partial^{2}I_{i}}{\partial x^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\cos^{2}\theta_{i}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\sin^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{\partial y^{\prime 2}} = {{\frac{\partial^{2}I_{i}}{\partial x^{2}}\sin^{2}\theta_{i}} - {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{\partial y^{2}}\cos^{2}\theta_{i}}}},{\frac{\partial^{2}I_{i}}{{\partial x^{\prime}}{\partial y^{\prime}}} = {{{- \frac{1}{2}}\frac{\partial^{2}I_{i}}{\partial x^{2}}{\sin \left( {2\theta_{i}} \right)}} + {\frac{\partial^{2}I_{i}}{{\partial x}{\partial y}}{\cos \left( {2\theta_{i}} \right)}} + {\frac{1}{2}\frac{\partial^{2}I_{i}}{\partial y^{2}}{\sin \left( {2\theta_{i}} \right)}}}}$